se <- function(x) {sd(x)/sqrt(length(x))}

gov.cohesion <- function(mat) {
    if(nrow(mat)==1) {
        return(0)
    } else {
        return(max(dist(mat)))
    }
}

utility.fun <- function(alpha, dim1, dim2, pprefdim1, pprefdim2, size, gov.seats, discount, govmember, caretaker) {
                         policy <- alpha*(1-((dim1-pprefdim1)^2)-((dim2-pprefdim2)^2))
                         office <- (1-alpha)*((govmember*size)/gov.seats)
                         return((policy+office)*(discount^caretaker))
                        }



eff.num.gov.parties <- function(size, gov.seats) {
    shares <- size/gov.seats
    return(1/sum(shares*shares))
}


vector.euclidean.distancer <- function(input) {
    sqrt(((input[1]-input[3])^2)+((input[2]-input[4])^2))
}



euclidean.distancer <- function(moving, target1, target2) {
    sqrt(((moving[1]-target1)^2)+((moving[2]-target2)^2))
}

find.gov.position <- function(dim1, dim2, seats, govmember) {
    sumseats <- sum(seats[govmember==1])
    gov1 <- sum(dim1[govmember==1]*seats[govmember==1])/sumseats
    gov2 <- sum(dim2[govmember==1]*seats[govmember==1])/sumseats
    return(c(gov1, gov2))
}

new.gov.represener <- function(vpref, gov1, gov2) {
    gmedian <- l1median(vpref)
    sum2.dist.voter.median.voter <- sum(apply(vpref, 1, euclidean.distancer,
                                          gmedian[1], gmedian[2])^2)
    distance2.to.gov <- apply(vpref, 1, euclidean.distancer, unique(govdim1), unique(govdim2))^2    
}


make.ensemble <- function(i) { #loop over all parts of combined results
    if (i > 5000) {
  ensemble <- readRDS(paste0("/beegfs/work/ws/ma_rlehrer-store-0/no_evo/results/sim_", i, ".rds"))
  govpos <- unique(ensemble[,c("gov1", "gov2"):=as.list(find.gov.position(dim1,
                                                  dim2, size, govmember)),
                            list(num, round)][,list(num, round, gov1, gov2)])
  
  voter <- readRDS(paste0("/beegfs/work/ws/ma_rlehrer-store-0/no_evo/voters/voters_", i, ".rds"))
  voter$experiment <- NULL
  voter <- voter[,c("l1med1", "l1med2"):=as.list(l1median(cbind(vprefdim1, vprefdim2))), num]
  voter$dist2.voter.median <- apply(voter[,list(vprefdim1, vprefdim2, l1med1, l1med2)], 1,
                                    vector.euclidean.distancer)^2
  sum.dist2.voter.median <- unique(voter[,sum.dist2.voter.median:=sum(dist2.voter.median), num][,num, sum.dist2.voter.median])
  voter <- voter[,list(num, vprefdim1, vprefdim2)]

  voter <- merge(govpos, voter, by=c("num"), all.y=TRUE, allow.cartesian=TRUE)
  voter$dist2.voter.gov <- sqrt(((voter$vprefdim1-voter$gov1)^2)+((voter$vprefdim2-voter$gov2)^2))^2
  voter <- unique(voter[,sum.dist2.voter.gov:=sum(dist2.voter.gov), list(num, round)][,list(num, round, sum.dist2.voter.gov)])

  rep <- merge(voter, sum.dist2.voter.median, by="num", all.x=TRUE, allow.cartesian=TRUE)
  rep$gov.representation2 <- 1-(rep$sum.dist2.voter.median/rep$sum.dist2.voter.gov)
  rep <- rep[,list(num, round, gov.representation2)]
  ensemble <- merge(ensemble, rep, by=c("num", "round"))
  rm(rep)
  rm(voter)
  rm(govpos)
  rm(sum.dist2.voter.median)

  ensemble <- ensemble[,gov.seats:=sum(size[govmember==1]),
                       list(num, experiment, round)]
  ensemble <- ensemble[,c("num.gov.parties", "gov.cohesion", "oversized", "eff.num.gov.parties"):=
                       list(sum(govmember),
                            gov.cohesion(cbind(dim1[govmember==1], dim2[govmember==1])),
                            max(ifelse((gov.seats[govmember==1]-size[govmember==1])>500,1,0)),
                            eff.num.gov.parties(size[govmember==1], gov.seats[govmember==1])
                            ),  list(num, experiment, round)]
  ensemble <- ensemble[,majority:=ifelse(gov.seats>500, 1, 0)]
  ensemble <- ensemble[,utility:=utility.fun(alpha, dim1, dim2, pprefdim1, pprefdim2, size, gov.seats, discount, govmember, caretaker)]
  ensemble <- ensemble[,mean.utility:=mean(utility), list(experiment, num, round)]
  
 #keep only one obs per system for system-level vars
 ensemble <- ensemble[party==1, list(experiment, round, rule, alpha, memory, aspiration,
                                parties, discount, sd.distance, representation,
                                eccentricity, gov.representation, gov.representation2,
                                gov.seats, num.gov.parties, gov.cohesion, oversized, majority,
                                eff.num.gov.parties, caretaker, mean.utility
                          )]
 
 #compute ensemble mean of system-vars
 ensemble <- unique(ensemble[,c("mean.eccentricity", "mean.representation",
                                  "mean.gov.representation", "mean.gov.representation2",
                                  "se.eccentricity", "se.representation",
                                  "se.gov.representation", "se.gov.representation2",
 "mean.gov.seats", "mean.num.gov.parties", "mean.gov.cohesion", "mean.oversized", "mean.majority",
 "se.gov.seats", "se.num.gov.parties", "se.gov.cohesion", "se.oversized", "se.majority",
 "mean.ENgov", "se.ENgov", "mean.caretaker", "se.caretaker", "mean.mean.utility", "sd.mean.utility") :=
              list(mean(eccentricity),
                   mean(representation),
                   mean(gov.representation),
                   mean(gov.representation2),                   
                    se(eccentricity),
                   se(representation),
                   se(gov.representation),
                   se(gov.representation2),                   
                   mean(gov.seats),
                   mean(num.gov.parties),
                   mean(gov.cohesion),
                   mean(oversized),
                   mean(majority),
                   se(gov.seats),
                   se(num.gov.parties),
                   se(gov.cohesion),
                   se(oversized),
                   se(majority),
                   mean(eff.num.gov.parties),
                   se(eff.num.gov.parties),
                   mean(caretaker),
                   se(caretaker),
                   mean(mean.utility),
                   sd(mean.utility)),
              list(experiment, round)][,list(experiment,
                                                        round,
                                                        rule,
                                                        alpha,
                                                        memory,
                                                        aspiration,
                                                        parties,
                                                        discount,
                                                        sd.distance,
                                 mean.eccentricity,
                                 mean.representation,
                                 mean.gov.representation,
                                 mean.gov.representation2,                                          
                                 se.eccentricity,
                                 se.representation,
                                 se.gov.representation,
                                 se.gov.representation2,                                            
                                 mean.gov.seats,
                                 mean.num.gov.parties,
                                 mean.gov.cohesion,
                                 mean.oversized,
                                 mean.majority,
                                 se.gov.seats,
                                 se.num.gov.parties,
                                 se.gov.cohesion,
                                 se.oversized,
                                 se.majority,
                                 mean.ENgov,
                                 se.ENgov,
                                 mean.caretaker,
                                 se.caretaker,
                                             mean.mean.utility,
                                             sd.mean.utility)])
 saveRDS(ensemble, file=paste0("/beegfs/work/ws/ma_rlehrer-store-0/no_evo/ensemble/ensemble",i,".rds"))
}

} else {
###########################
# HUTER ENSEMBLE AVERAGES #
###########################
  ensemble <- readRDS(paste0("/beegfs/work/ws/ma_rlehrer-store-0/no_evo/results/sim_", i, ".rds"))
  govpos <- unique(ensemble[,c("gov1", "gov2"):=as.list(find.gov.position(dim1,
                                                  dim2, size, govmember)),
                            list(num, round)][,list(num, round, gov1, gov2)])
  
  voter <- readRDS(paste0("/beegfs/work/ws/ma_rlehrer-store-0/no_evo/voters/voters_", i, ".rds"))
  voter$experiment <- NULL
  voter <- voter[,c("l1med1", "l1med2"):=as.list(l1median(cbind(vprefdim1, vprefdim2))), num]
  voter$dist2.voter.median <- apply(voter[,list(vprefdim1, vprefdim2, l1med1, l1med2)], 1,
                                    vector.euclidean.distancer)^2
  sum.dist2.voter.median <- unique(voter[,sum.dist2.voter.median:=sum(dist2.voter.median), num][,num, sum.dist2.voter.median])
  voter <- voter[,list(num, vprefdim1, vprefdim2)]

  voter <- merge(govpos, voter, by=c("num"), all.y=TRUE, allow.cartesian=TRUE)
  voter$dist2.voter.gov <- sqrt(((voter$vprefdim1-voter$gov1)^2)+((voter$vprefdim2-voter$gov2)^2))^2
  voter <- unique(voter[,sum.dist2.voter.gov:=sum(dist2.voter.gov), list(num, round)][,list(num, round, sum.dist2.voter.gov)])

  rep <- merge(voter, sum.dist2.voter.median, by="num", all.x=TRUE, allow.cartesian=TRUE)
  rep$gov.representation2 <- 1-(rep$sum.dist2.voter.median/rep$sum.dist2.voter.gov)
  rep <- rep[,list(num, round, gov.representation2)]
  ensemble <- merge(ensemble, rep, by=c("num", "round"))
  rm(rep)
  rm(voter)
  rm(govpos)
  rm(sum.dist2.voter.median)
 ensemble <- ensemble[,gov.seats:=sum(size[govmember==1]),
                       list(experiment,round)]

  ensemble <- ensemble[,c("num.gov.parties", "gov.cohesion", "oversized", "eff.num.gov.parties"):=
                       list(sum(govmember),
                            gov.cohesion(cbind(dim1[govmember==1], dim2[govmember==1])),
                            max(ifelse((gov.seats[govmember==1]-size[govmember==1])>500,1,0)),
                            eff.num.gov.parties(size[govmember==1], gov.seats[govmember==1])
                            ),  list(experiment,round)]
  ensemble <- ensemble[,majority:=ifelse(gov.seats>500, 1, 0)]
  ensemble <- ensemble[,utility:=utility.fun(alpha, dim1, dim2, pprefdim1, pprefdim2, size, gov.seats, discount, govmember, caretaker)]
  ensemble <- ensemble[,mean.utility:=mean(utility), list(experiment, num, round)]


 #keep only one obs per system for system-level vars
 ensemble <- ensemble[party==1, list(experiment, round, rule, alpha, memory, aspiration,
                                parties, discount, sd.distance, representation,
                                eccentricity, gov.representation, gov.representation2,
                                gov.seats, num.gov.parties, gov.cohesion, oversized, majority,
                                eff.num.gov.parties, caretaker, mean.utility
                          )]


 #compute ensemble mean of system-vars
 ensemble <- unique(ensemble[round>200,c("mean.eccentricity", "mean.representation",
                                  "mean.gov.representation", "mean.gov.representation2",
                                  "se.eccentricity", "se.representation",
                                  "se.gov.representation", "se.gov.representation2",
 "mean.gov.seats", "mean.num.gov.parties", "mean.gov.cohesion", "mean.oversized", "mean.majority",
 "se.gov.seats", "se.num.gov.parties", "se.gov.cohesion", "se.oversized", "se.majority",
 "mean.ENgov", "se.ENgov", "mean.caretaker", "se.caretaker", "mean.mean.utility", "sd.mean.utility") :=
              list(mean(eccentricity),
                   mean(representation),
                   mean(gov.representation),
                   mean(gov.representation2),                   
                    se(eccentricity),
                   se(representation),
                   se(gov.representation),
                   se(gov.representation2),
                   mean(gov.seats),
                   mean(num.gov.parties),
                   mean(gov.cohesion),
                   mean(oversized),
                   mean(majority),
                   se(gov.seats),
                   se(num.gov.parties),
                   se(gov.cohesion),
                   se(oversized),
                   se(majority),
                   mean(eff.num.gov.parties),
                   se(eff.num.gov.parties),
                   mean(caretaker),
                   se(caretaker),
                   mean(mean.utility),
                   sd(mean.utility)),
              list(experiment)][,list(experiment,
                                                        round,
                                                        rule,
                                                        alpha,
                                                        memory,
                                                        aspiration,
                                                        parties,
                                                        discount,
                                                        sd.distance,
                                 mean.eccentricity,
                                 mean.representation,
                                 mean.gov.representation,
                                 mean.gov.representation2,                                      
                                 se.eccentricity,
                                 se.representation,
                                 se.gov.representation,
                                 se.gov.representation2,                                      
                                 mean.gov.seats,
                                 mean.num.gov.parties,
                                 mean.gov.cohesion,
                                 mean.oversized,
                                 mean.majority,
                                 se.gov.seats,
                                 se.num.gov.parties,
                                 se.gov.cohesion,
                                 se.oversized,
                                 se.majority,
                                 mean.ENgov,
                                 se.ENgov,
                                 mean.caretaker,
                                 se.caretaker,
                                      mean.mean.utility,
                                      sd.mean.utility)])

ensemble <- unique(ensemble[round==601,
                            list(experiment,
                                 round,
                                 rule,
                                 alpha,
                                 memory,
                                 aspiration,
                                 parties,
                                 discount,
                                 sd.distance,
                                 mean.eccentricity,
                                 mean.representation,
                                 mean.gov.representation,
                                 mean.gov.representation2,                                 
                                 se.eccentricity,
                                 se.representation,
                                 se.gov.representation,
                                 se.gov.representation2,
                                 mean.gov.seats,
                                 mean.num.gov.parties,
                                 mean.gov.cohesion,
                                 mean.oversized,
                                 mean.majority,
                                 mean.ENgov,
                                 se.gov.seats,
                                 se.num.gov.parties,
                                 se.gov.cohesion,
                                 se.oversized,
                                 se.majority,
                                 se.ENgov,
                                 mean.caretaker,
                                 se.caretaker,
                                 mean.mean.utility,
                                      sd.mean.utility)])
ensemble$round <- 76
saveRDS(ensemble, file=paste0("/beegfs/work/ws/ma_rlehrer-store-0/no_evo/ensemble/ensemble", i, ".rds")
} #end else
} #end function

## ensemble <- list()
## for (i in c(1, 51:250)) {
##  ensemble[[i]] <- readRDS(paste0("/beegfs/work/ws/ma_rlehrer-mc_newfor_no_evo_results-0/ensemble/ensemble",i,".rds"))
## }
## ensemble <- rbindlist(ensemble)
## ## setnames(ensemble, "representation.mean", "mean.representation")
## ## setnames(ensemble, "eccentricity.mean", "mean.eccentricity")
## ## setnames(ensemble, "gov.representation.mean", "mean.gov.representation")
## ## setnames(ensemble, "representation.se", "se.representation")
## ## setnames(ensemble, "eccentricity.se", "se.eccentricity")
## ## setnames(ensemble, "gov.representation.se", "se.gov.representation")
## ## ensemble <- melt(ensemble, measure=patterns("^mean.", "^se."), value.name=c("mean", "se"))
## ## ensemble$variable <- factor(ensemble$variable,
## ##                             labels=c("eccentricity.mean","representation.mean", "gov.rep.mean"))
## save(ensemble, file="/beegfs/work/ws/ma_rlehrer-mc_newfor_no_evo_results-0/analysis/ensemble.Rdata")

